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In the large-[7 limit, the Falicov-Kimball model maps onto an effective Ising model, with an order 
parameter described by a BCS-like mean-field theory in infinite dimensions. In the small-f7 limit, 
van Dongen and Vollhardt showed that the order parameter assumes a strange non-BCS-like shape 
with a sharp reduction near T « T c /2. Here we numerically investigate the crossover between these 
two regimes and qualitatively determine the order parameter for a variety of different values of 
U . We find the overall behavior of the order parameter as a function of temperature to be quite 
anomalous. 



PACS numbers: 71.10.-W, 71.30.+h, 71.45.-d, 72.10.-d 

Dynamical mean field theoryi has been widely used to 
study electron correlations in a variety of different inter- 
acting systems. Much work has focused on the "param- 
agnetic" metal-insulator transition, on transport proper- 
ties in the normal state, and on determining phase dia- 
grams to ordered phases via a susceptibility analysis or 
a Maxwell construction. The properties of the ordered 
phase have been less studied, yet there is much interesting 
physics to examine there. For example, one might have 
thought that since the system is infinite-dimensional, 
both the critical behavior and the order parameter as a 
function of temperature would be determined by a BCS- 
like mean-field picture. Indeed, the critical exponents are 
always mean-field like, but the order parameter can have 
quite anomalous behavior as a function of temperature. 
This anomalous behavior is amplified for small correla- 
tion strength, since many models map onto effective spin 
models for large correlations, and the order parameter of 
an infinite-dimensional spin model is always mean-field- 
like. In this contribution, we examine in detail the case 
of a commensurate (two-sublattice) charge-density-wave 
(CDW) phase in the spinless Falicov-Kimball model at 
half filling. Analytical work in the large and small- J7 
limits has already been carried out2*^. 

The Falicov Kimball (FK) modeU was introduced in 
1969 to describe metal-insulator transitions in transition- 
metal and rare-earth compounds. The model consists 
of two types of particles: itinerant conduction electrons 
and localized ions that mutually interact with an on-site 
Coulomb interaction. It is the simplest Fermionic model 
of electron correlations and the spinless version can be 
interpreted as a simple model for crystallization, 5 where 
the system has a phase transition from a disordered phase 
at high temperature to an ordered phase as the temper- 
ature is lowered. We examine the order parameter as a 
function of temperature here. 

Many body effects enter via the statistical mechanics 
associated with annealed averaging. The FK model is the 
simplest many-body problem that can be solved exactly 



in the limit of large dimensions. Brandt and MielscbJliLS 
presented the first solution of this problem using dynam- 
ical mean-field theory. Their solution quantitatively il- 
lustrated how a period-two CDW phase is stabilized at 
low temperatures. 

For the symmetric half-filled case on a bipartite lat- 
tice, the ordering is into a commensurate CDW state, in 
which the particles order in a checkerboard pattern: the 
conduction electrons preferentially occupy one sublattice, 
and the ions the other one. Van Dongen and Vollhardt2^ 
derived an analytical expression for the critical temper- 
ature T c as a function of the Coulomb interaction in the 
large and small- U limits. For large- U, it has the conven- 
tional strong-coupling t 2 /U form, but for small coupling, 
T c is much larger than the BCS-like exponential form of 
exp[— C/U\. They obtained similar results for the order 
parameter — for large-£/ it behaved mean-field-like, but 
for small-f7 the behavior was quite anomalous with a 
sharp decrease as T — > T c /2. Here we numerically ex- 
plore the behavior of the order parameter between these 
two limits. A self-consistent algorithm^ is employed to 
perform efficient computations at low temperature and to 
make it feasible to study the order parameter when U is 
small. Surprisingly, the order parameter is not BCS-like! 

The spinless Falicov-Kimball model is represented by 
the following Hamiltonian: 



Ef^flfi + u^ctartfi. (l) 



The conduction electrons (created or destroyed at site i 
by c| or Cj) can hop between nearest-neighbor sites with 



a hopping matrix 1 —tij = —t*/2yD. The localized ions 
(created or destroyed at site i by fj or /,*) have a site 
energy Ef. There is a Coulomb interaction U between 
the localized ions and the conduction electrons that sit at 
the same lattice site. A chemical potential \i is employed 
to determine the number of conduction electrons. At half 
filling, we have fi = U/2 and Ef — —U/2. 
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We examine two bipartite lattices here — the hypercu- 
bic lattice with a density of states (DOS) given by p(e) = 
exp(— e 2 )/t* v / 7r and the infinite-coordination Bethe lat- 
tice with p(e) = V4 — e 2 /t*2ir (we take t* — 1 as our 
energy unit). The local Green's function G{iuj n ) can be 
written as the Hilbert transform of the noninteracting 
DOSp(e) 



G(iu> n ) = / dep(e) 



1 



H - S, 



(2) 



where u> n = nT(2n + 1) is the Fermionic Matsubara fre- 
quency and Yj(iu) n ) = S„ is the local self energy. Dyson's 
equation for the local self energy reads 



iuj n + p - A„ - G 1 (iuj n ) 



(3) 



with A n the dynamical mean field evaluated at the nth 
Matsubara frequency. The effective medium Green's 
function Go(ico n ) is 



G (iu n ) 



1 



1 



G 1 (iw n ) + E„ ioj n + p,- X, 



(4) 



Solving the atomic problem in a time-dependent field 
yields another equation for the local Green's function 



G(iuj n ) 



1 — W-y 



!"1 



iw n + p- A„ iuj n + p - A„ - U 



(5) 



where w± represents the density of the ions which is de- 
termined by the atomic partition function 



Z at {\) = Z (X,p) + e~ l3E fZ (X,p-U) (6) 



with 



Z (A,/i) = 2e^ J] 



(7) 



and wi = e-P E fZo{X,p- U)/Z at . When E f = -U/2 
and p = U/2, then A„ is pure imaginary and we have 
wi=0.5. 

The iterative algorithm for finding the self-energy is 
that of JarrelUS: (i) first set the self-energy equal to 
zero; (ii) then determine the local Green's function from 
Eq. (J2J ; (iii) then determine the effective-medium Green's 
function from Eq. 10} and the ion density from the atomic 
partition function (in our case we always have w% = 0.5); 
(iv) next determine the new local Green's function from 
Eq. (JSJ and the new self energy from Eq. J3J). Starting 
with the new self-energy, repeat steps (ii-iv) until con- 
vergence is reached. 

After solving for the self-energy and the Green's func- 
tions, the susceptibility to a CDW instability (T\ = 
{niUj) — (ni)(rij}) can be calculated following the deriva- 
tion of Brandt and Mielsch^i At high temperature, the 
electrons are uniformly distributed throughout the lat- 
tice. As the temperature is lowered, the charge density 
becomes nonuniform when the momentum-dependent 



susceptibility diverges and a CDW forms. After some 
tedious algebra, the transition temperature is found to 
occur when 



1 



^ (l + G n E n )[l + G„(E n -E0] 

1 



1 + G„[E„ - (1 - Wl )U] + G nVn (q)[l + G„(2E, 
holds, where 

G r , 



G„ 1 



(9) 



and the bare susceptibility satisfies 

in terms of the momentum-dependent Green's functions 
(N is the number of lattice sites). The q-dependence of 
X can be summarized by the scalar parameter 11 



cos qi . 



(11) 



For X = — 1, there is a phase transition from a uniform 
distribution to a (two-sublattice) checkerboard phase 
when the temperature falls below T c . This transition 
temperature is plotted versus U for (a) the Bethe lattice 
and (b) the hypercubic lattice in Fig. 




0.2 0.4 0.6 0.8 

Interaction Strength U/(t*+U) 

FIG. 1: Transition temperature for the (X — —1) CDW at 
half filling on (a) the Bethe lattice and (b) the hypercubic 
lattice as determined from the solution of Eq. (0 . The hori- 
zontal axis is plotted as U/(t* + U) in order to include both 
the weak and strong coupling results on the same graph. 
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For T < T c , we need to revise the algorithm to solve for 
the Green's functions and to determine the order param- 
eter. Since the thermodynamic quantities now differ on 
the A and B sublattices of the bipartite lattice, we label 
them by their respective sublattice. The order parameter 
for the CDW is defined to be A = \wf - wf |. We find 
that if we try to iterate a generalization of the Jarrell 
algorithm to include the possibility of CDW order, then 
the equations suffer from "critical slowing down" as one 
approaches T c from below. This would make an efficient 
determination of the order parameter as a function of T 
impossible when T c became small because the number 
of Matsubara frequencies employed becomes large, and 
the equations must be iterated for a long time to achieve 
convergence. Instead, we adopt a variant of the algo- 
rithm proposed by Gruber et al&. We start by setting 
the self energies on each sublattice to zero. Then we (i) 
choose a value for the order parameter < A < 1 and set 
= 0.5(1 + A) and wf = 0.5(1 - A); (ii) we determine 
the parameter Z n defined by 

Zn = ^Ki/i-S^^ + ^-SS) (12 ) 

and calculate the local Green's function on the A sublat- 
tice 



finally (ix) we calculate the new self energy 



Gt = Wn + t! S " f de- P{e) 



(13) 



Z n , J Z n — e 

(iii) we solve for the effective medium on the A sublattice 

Gfoiu, n ) = (G£- 1 + V*)- 1 ; (14) 
and (iv) we calculate the new local Green's function 

Gi = (1 - wt)G$(iu n ) + — j^- ; (15) 

Gq (iuj n ) - U 



finally (v) we determine the new self energy 



(16) 



Next we perform a similar analysis on the B sublattice: 
(vi) first we update Z n then calculate the local Green's 
function on the i?-sublattice 



u _ itjJ n + fi - f p{e) 



G„ — 



(17) 



Z n J Z„ — e 

(vii) we determine the effective medium 

GZ(iu n ) = (G2- 1 +Z*)- 1 ; (18) 
and (viii) we find the new local Green's function 



G» = (l-OG£ («"„) + 



Gq (iuj n )-U 



(19) 



l , . 



(20) 



Steps (ii-ix) are repeated until converged. Then we (x) 
calculate the local ion site energy on each sublattice by 
solving the generalization of the formula for the filling 
wi to the two-sublattice case, and solving for the Ef on 
each sublattice: 



Ef = 



U 



and 



Ef 



U 



-Tin 



-Tin 



1 



-Tj2W-UG^(ico n )}, (21) 



-Tj2Hl-UG$(iuj n )]. (22) 



Ej , then a consistent thermodynamic solution 



If Ef 

has been achieved. If not, then we adjust the order pa- 
rameter A and repeat the algorithm starting at step (i) 
until a consistent solution is reached. Since the ion den- 
sities on each sublattice are not updated during the iter- 
ative part of the algorithm, this technique does not suffer 
from critical slowing down, and it is easy to perform on 
both the Bethe lattice and the hypercubic lattice, sim- 
ply by changing the corresponding DOS in the Hilbcrt 
transforms of Eqs. I1K'>[1 and i|17|) . In our numerical cal- 
culations, we use an energy cutoff of lOt* for the Mat- 
subara frequencies; the maximum number of Matsubara 
frequencies we include is 250 000. 

For U ~ * oo, the FK model maps onto an antiferromag- 
netic Ising model with J = t* 2 /4DU. This spin system 
can be analyzed with static mean field theory to yield 
the order parameter in a Curie- Weiss form, 



AT 

A = tanh(^), 



(23) 



with T c = t* 2 /4U on the hypercubic lattice. For U — ► 0, 
the order parameter varies sharply from this result, with 
a steep decrease near T c /2 and a flattening out before 
going to zero at T c . We examine the transition between 
these two limits numerically. We find that U must be 
very small to even see the deviations from the mean- 
field theory curve. The evolution from the small-?/ to 
large-£7 limits is nonmonotonic. As U increases, the order 
parameter curve becomes steeper and steeper near T c 
and flatter and flatter for small T until we reach the U 
corresponding to approximately the maximum of the T c 
versus U curve. Once U increases past the maximum, 
then the order parameter curve slowly approaches the 
large- U limit. 

We illustrate this behavior in Fig.|21for the Bethe lat- 
tice. When U lies below 0.2 [panel (a)], we can see the 
depression in the curve develop as we evolve toward the 
U = limit. For 0.2 < U < 2.6, the curve moves the 
right and becomes steeper (but the critical exponent al- 
ways remains 0.5) [panel (b)]. Finally, for larger values 
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0.2 0.4 0.6 0.8 1 

FIG. 2: Order parameter for the (X = -1) CDW of the FK 
model at half filling on the Bethe lattice. Panel (a) shows 
the small-(7 results where a depression can be seen (U — * 0, 
0.005, 0.03, 0.1, and 0.25), panel (b) shows the intermediate 
results ((7=0.25, 0.5, and 2.6), where the order parameter 
curve moves out to the right as U is increased, and panel (c) 
shows the large-?/ regime ((7=2.6, 3.0, 4.0 and oo), where the 
curve shifts to the left and ultimately takes the mean-field- 
theory form. 



body systems, to the infinite-dimensional limit, or to the 
FK model. We tend to believe that it is most likely a cu- 
rious property of the FK model itself. The weak-coupling 
limit of the FK model is known to have anomalous be- 
havior for it's T c , being much higher than what would be 
predicted by a BCS approach, and it is possible that the 
anomalous shape survives in finite dimensions as well. 




0.2 0.4 0.6 0.8 1 



FIG. 3: Order parameter for the (X = -1) CDW of the 
FK model at half filling on the hypercubic lattice. Panel (a) 
shows the small-!7 results where a depression can be seen 
((7= 0.005, 0.03, 0.1, 0.5, and 1.0) and panel (b) shows the 
larger-(7 results ((7=1.0, 2.0, 4.0, and oo), where the order 
parameter curve moves back to the left as U is increased, but 
overshoots, and then moves to the right to ultimately take 
the mean-field-theory form. 



of U, the curve moves to the left until it assumes the 
mean- field-theory form [panel (c)]. 

The results on the hypercubic lattice are similar and 
are shown in Fig. El We show only two different regimes, 
as the curve is steepest when U lies near 1.0. Note that 
there are no analytic results available for the — > limit 
on the hypercubic lattice. The one difference from the 
Bethe lattice, is that here, the curves move to the right 
as we increase from U = 1 to U = 4, but they overshoot 
the mean-field theory curve, so they move back to the 
left to reach the U — ► oo limit for larger values of U. 

In summary, we have shown that the order parameter 
has an anomalous evolution as a function of tempera- 
ture for small U in the spinless FK model. One can ask 
is this anomalous behavior something generic to many- 



Unfortunately, our numerical calculations show that U 
must be very small in order to see the anomalous shape 
emerge, and it is likely that the T c is so small, that nu- 
merical calculations in finite dimensions may have trou- 
ble being performed at the temperatures necessary to see 
the anomaly. Even in infinite-dimensions, we are limited 
by how small we can make U because the T c is pushed so 
low that the number of Matsubara frequencies included 
in the calculation grows too large. 
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